Projekt Seri Kohore

1. Hyrje

library(readxl)
data <- read_excel("Dataset3.xlsx")
head(data)

## # A tibble: 6 × 4
##    Viti ProdhimiDritherave PlehraKimike   CO2
##   <dbl>              <dbl>        <dbl> <dbl>
## 1  1961             293932         6650  1.29
## 2  1962             313508         4110  1.35
## 3  1963             293449         4040  1.10
## 4  1964             350144         3990  1.03
## 5  1965             335293         3830  1.08
## 6  1966             433094         5560  1.23

Qëllimi i Projektit:

Ky projekt synon të analizojë ndryshimet dhe trendet e prodhimit të drithërave në lidhje me përdorimin e plehrave kimikë dhe emëtimet e CO2 në një periudhë afatgjatë, duke përdorur analiza të serive kohore.

Pyetjet Kërkimore:

Si ka ndryshuar prodhimi i drithërave gjatë periudhës së studimit?
A ka një ndikim të drejtpërdrejtë përdorimi i plehrave kimikë dhe emetimet e CO2 në prodhimin e drithërave?
Si mund të parashikohen ndryshimet e ardhshme në prodhimin e drithërave?

2. Përgatitja e të Dhënave

Përshkrimi i Të Dhënave:

Variablat: - Viti: Viti i të dhënave. - ProdhimiDritherave: Prodhimi i drithërave (kg për hektar). - PlehraKimike: Përdorimi i plehrave kimikë (kg për hektar). - CO2: Emetimet e CO2 .

Burimi i të Dhënave:

Këto të dhëna janë të besueshme dhe janë burim nga “Our World in Data” ,të cilat përfaqësojnë një periudhë 59-vjeçare (1961 - 2020), duke lejuar analizën e tendencave afatgjata dhe të ndikimeve të ndërsjella midis variablave të përzgjedhur.

Konvertimi i datave në format kohor:

data$Date <- as.Date(paste0(data$Viti, "-01-01"))

Kontrolli për vlera që mungojnë:

missing_values <- sum(is.na(data))
cat("Numri i vlerave që mungojnë në dataset është:", missing_values, "\n")

## Numri i vlerave që mungojnë në dataset është: 0

Meqë numri i vlerave të munguara është 0, do të thotë që dataset-i është i plotë dhe nuk kemi nevojë të merremi me trajtimin e mungesave në të dhëna.Tani mund të kalojmë direkt në analizat e tjera.

3. Analiza Eksploruese e të Dhënave

3.1 Vizualizimet e serive : 3.1.1) Vizualizimi i serisë kohore për Prodhimin e Drithërave

library(ggplot2)
ggplot(data, aes(x = Date, y = ProdhimiDritherave)) +
  geom_line(color = "blue", size = 1) +
  ggtitle("Prodhimi i Drithërave nga Viti 1961 në 2020") +
  scale_x_date(date_labels = "%Y", date_breaks = "4 years") +
  labs(x = "Viti", y = "Prodhimi i Drithërave") +
  theme_minimal() +
  theme(
    axis.text.x = element_text(angle = 45, hjust = 1),
    plot.title = element_text(hjust = 0.5, size = 14, face = "bold")
  )

## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

Grafiku i serisë kohore tregon një trend rritës afatgjatë, me disa luhatje gjatë viteve 1960-1980. Pas vitit 1990, rritja bëhet më e qëndrueshme.

3.1.2) Vizualizimi i serisë kohore për CO2

ggplot(data, aes(x = Date, y = CO2)) +
  geom_line(color = "green", size = 1) +
  ggtitle("Emëtimi i CO2 nga Viti 1961 në 2020") +
  scale_x_date(date_labels = "%Y", date_breaks = "4 years") +
  labs(x = "Viti", y = "CO2 (metric tonnes)") +
  theme_minimal() +
  theme(
    axis.text.x = element_text(angle = 45, hjust = 1),
    plot.title = element_text(hjust = 0.5, size = 14, face = "bold")
  )

CO2 ka një trend rritës të qëndrueshëm, që mund të lidhet me rritjen e industrializimit dhe aktivitetet intensive bujqësore.Pas vitit 1995, rritja bëhet më e qëndrueshme.

3.1.3) Vizualizimi i serisë kohore për Plehrat Kimike

ggplot(data, aes(x = Date, y = PlehraKimike)) +
  geom_line(color = "red", size = 1) +
  ggtitle("Plehrat Kimike nga Viti 1961 në 2020") +
  scale_x_date(date_labels = "%Y", date_breaks = "4 years") +
  labs(x = "Viti", y = "Plehra Kimike (tonë)") +
  theme_minimal() +
  theme(
    axis.text.x = element_text(angle = 45, hjust = 1),
    plot.title = element_text(hjust = 0.5, size = 14, face = "bold")
  )

Grafiku për plehrat kimike ka një rritje të fortë lineare nga vitet 1961 deri në 2020, që tregon një intensifikim të përdorimit të plehrave kimike për të përmirësuar produktivitetin bujqësor.Pas vitit 1990, rritja bëhet më e qëndrueshme.

3.2 Statistika përmbledhëse:

summary(data)

##       Viti      ProdhimiDritherave  PlehraKimike         CO2        
##  Min.   :1961   Min.   : 293449    Min.   :  3830   Min.   :0.4742  
##  1st Qu.:1976   1st Qu.: 503410    1st Qu.: 21522   1st Qu.:1.2495  
##  Median :1990   Median : 645713    Median : 49485   Median :1.6217  
##  Mean   :1990   Mean   : 660907    Mean   : 50074   Mean   :1.5921  
##  3rd Qu.:2005   3rd Qu.: 741402    3rd Qu.: 62268   3rd Qu.:1.8938  
##  Max.   :2020   Max.   :1067000    Max.   :113010   Max.   :2.7379  
##       Date           
##  Min.   :1961-01-01  
##  1st Qu.:1975-10-01  
##  Median :1990-07-02  
##  Mean   :1990-07-02  
##  3rd Qu.:2005-04-02  
##  Max.   :2020-01-01

3.2.1) Vlera mesatare për seritë

mesataret <- sapply(data[, c("ProdhimiDritherave", "PlehraKimike", "CO2")], mean)
print(mesataret)

## ProdhimiDritherave       PlehraKimike                CO2 
##       6.609074e+05       5.007350e+04       1.592133e+00

3.2.2) Devijimi standart për seritë

devijimistandart <- sapply(data[, c("ProdhimiDritherave", "PlehraKimike", "CO2")], sd)
print(devijimistandart)

## ProdhimiDritherave       PlehraKimike                CO2 
##       2.095177e+05       3.314175e+04       5.783174e-01

Variablat tregojnë shpërndarje të gjerë, sidomos për prodhimin e drithërave dhe plehrave kimike.

3.3 Autokorrelacionet:

# Autokorrelacioni për Prodhimi i Drithërave
acf(data$ProdhimiDritherave, main = "Autokorrelacioni për Prodhimi i Drithërave")

# Autokorrelacioni për Plehra Kimike
acf(data$PlehraKimike, main = "Autokorrelacioni për Plehra Kimike")

# Autokorrelacioni për CO2
acf(data$CO2, main = "Autokorrelacioni për CO2")

Grafiku ACF për Prodhimin e Drithërave tregon: • Vlera të larta të autokorrelacionit për laget e para, duke sugjeruar një lidhje të fortë mes vlerave të njëpasnjëshme të serisë kohore. • Gradualisht, autokorrelacionet zvogëlohen, por disa vlera janë ende mbi kufirin dhe për laget e mëdha, që sugjeron prani të një trendi afatgjatë.

Interpretimi: • Këto rezultate janë tipike për seri kohore jo stacionare, siç është konfirmuar edhe nga testi i stacionaritetit (ADF). • Trendet afatgjata dhe struktura e lidhur me kohën sugjerojnë nevojën për diferencim për të eliminuar ndikimin e trendit.

Grafiku ACF për Plehrat Kimike tregon: • Një model të ngjashëm me prodhimin e drithërave, ku autokorrelacionet fillestare janë shumë të larta dhe pastaj gradualisht zvogëlohen. • Kjo tregon se plehrat kimike kanë një strukturë të fortë trendi, që ndikon te vlerat në periudha të njëpasnjëshme.

Interpretimi: • Këto rezultate tregojnë që përdorimi i plehrave kimike ndjek një model të qëndrueshëm të rritjes afatgjatë. • Diferencimi ose përdorimi i modeleve të trendit si Holt-Winters është i përshtatshëm për modelimin e kësaj serie.

Grafiku ACF për Emëtimet e CO2: • Tregon gjithashtu autokorrelacione të larta për laget fillestare, por me një rënie më të shpejtë krahasuar me dy seritë e tjera. • Kjo sugjeron një trend më të moderuar afatgjatë në krahasim me prodhimin e drithërave dhe plehrat kimike.

Interpretimi: • Emetimët e CO2 kanë një lidhje më të dobët afatgjatë, por ende ndjekin një strukturë trendi. • Ky variabël është më pak i lidhur me ciklet kohore të ngjashme me plehrat kimike ose prodhimin e drithërave.

3.4 Matrica e Korrelacionit:

library(corrplot)

## corrplot 0.95 loaded

cor_matrix <- cor(data[, c("ProdhimiDritherave", "PlehraKimike", "CO2")])
print(cor_matrix)

##                    ProdhimiDritherave PlehraKimike       CO2
## ProdhimiDritherave          1.0000000    0.8933734 0.7527944
## PlehraKimike                0.8933734    1.0000000 0.8596405
## CO2                         0.7527944    0.8596405 1.0000000

# Vizualizimi i matricës së korrelacionit
corrplot(cor_matrix, method = "circle", type = "upper", tl.col = "black", tl.srt = 45)

Korrelacion i lartë pozitiv midis plehrave kimike dhe prodhimit të drithërave (r = 0.89). Korrelacion i moderuar midis CO2 dhe prodhimit të drithërave (r = 0.75).

3.5 Paraqitja e të dhënave si objekt i serive kohore:

# Krijimi i serisë kohore
ts_prodhimi_dritherave <- ts(data$ProdhimiDritherave, start = c(1961), frequency = 1)
ts_plhera_kimike <- ts(data$PlehraKimike, start = c(1961), frequency = 1)
ts_CO2 <- ts(data$CO2, start = c(1961), frequency = 1)

3.6 Kontrolli i Stacionaritetit për serinë e Prodhimit të Drithërave: ADF test

library(tseries)

## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

adf.test(ts_prodhimi_dritherave) #p-value (0.5701) me e madhe se 0.05  (Seria jo stacionare)

## 
##  Augmented Dickey-Fuller Test
## 
## data:  ts_prodhimi_dritherave
## Dickey-Fuller = -2.0115, Lag order = 3, p-value = 0.5701
## alternative hypothesis: stationary

Seria e prodhimit të drithërave është jo stacionare siç e pamë dhe nga ndërtimi i autokorrelacioneve, prandaj do duhet të bejmë diferencimin e serisë për të eleminuar trendin.

3.7 Diferencimi i serisë të Prodhimit të Drithërave dhe kontrolli i stacionaritetit:

ts_diff <- diff(ts_prodhimi_dritherave) # Për stacionarizim
adf.test(ts_diff)

## 
##  Augmented Dickey-Fuller Test
## 
## data:  ts_diff
## Dickey-Fuller = -3.8399, Lag order = 3, p-value = 0.02277
## alternative hypothesis: stationary

kpss.test(ts_diff)

## Warning in kpss.test(ts_diff): p-value greater than printed p-value

## 
##  KPSS Test for Level Stationarity
## 
## data:  ts_diff
## KPSS Level = 0.20209, Truncation lag parameter = 3, p-value = 0.1

Pas diferencimit, seria bëhet stacionare (p < 0.05), duke sugjeruar një trend rritës që duhet larguar për modelim të suksesshëm.

3.8 Modeli i Regresit Linear për prodhimin e drithërave:

library(forecast)
model_prod_dritherave <- lm(ProdhimiDritherave ~ PlehraKimike , data = data)
summary(model_prod_dritherave)

## 
## Call:
## lm(formula = ProdhimiDritherave ~ PlehraKimike, data = data)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -191590  -65348    8489   72648  190053 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  3.781e+05  2.234e+04   16.93   <2e-16 ***
## PlehraKimike 5.648e+00  3.730e-01   15.14   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 94950 on 58 degrees of freedom
## Multiple R-squared:  0.7981, Adjusted R-squared:  0.7946 
## F-statistic: 229.3 on 1 and 58 DF,  p-value: < 2.2e-16

# Vizualizimi i regresionit
ggplot(data, aes(x = PlehraKimike, y = ProdhimiDritherave)) +
  geom_point() +
  geom_smooth(method = "lm", formula = y ~ x , se = TRUE, color = "blue") +
  labs(title = "Ndikimi i Plehrave Kimike në Prodhimin e Drithërave",
       x = "Plehra Kimike", 
       y = "Prodhimi i Drithërave") +
  theme_minimal()

# Saktësia
accuracy(model_prod_dritherave)

##                      ME     RMSE      MAE       MPE     MAPE      MASE
## Training set 1.1399e-11 93351.64 77313.57 -2.613464 13.65529 0.4670403

Modeli midis prodhimit të drithërave dhe plehrave kimike ka një R-squared = 0.79, duke treguar se 79% e variacionit në prodhim shpjegohet nga plehrat kimike.Prandaj ky model është i mirë për të parashikuar prodhimin e drithërave duke përdorur ndryshoren e plehrave kimike.

Rezultatet: - Plehrat kimike kanë ndikim të rëndësishëm pozitiv (p < 0.001). - Ky model është i thjeshtë, por i fuqishëm për interpretim.

4. Ndërtimi i Modelimi të Serive Kohore

4.1 Mesatarja e lëvizshme për prodhimin e drithërave:

library(TTR)
# Moving Average (MA) me periudhë 3
ma_data <- ma(ts_prodhimi_dritherave, order = 3)

plot(ts_prodhimi_dritherave, main = "Seri Kohore dhe Moving Average", col = "blue", type = "l")
lines(ma_data, col = "green", lwd = 2)
legend("topright", legend = c("Seri Kohore", "Moving Average"), col = c("blue", "green"), lwd = 2)

Ky grafik paraqet serinë kohore të prodhimit të drithërave (1961-2020) dhe mesataren e lëvizshme të saj me një periudhë prej 3 vitesh.
Vijat blu përfaqësojnë serinë kohore origjinale, ndërsa vija jeshile paraqet mesataren e lëvizshme.
Vijat blu tregojnë luhatjet e forta në prodhimin e drithërave, duke përfshirë disa luhatje të papritura, sidomos për vitet e hershme.
Mesatarja e lëvizshme (vija jeshile) zbut këto luhatje, duke dhënë një trend të qartë rritës të përgjithshëm gjatë periudhës.
Përdorimi i mesatares së lëvizshme është i dobishëm për të identifikuar trendin afatgjatë duke eliminuar luhatjet e rastësishme.
Shohim që vijat përshtaten mirë me njëra-tjetrën.
4.2 Mesatarja e lëvizshme eksponenciale për prodhimin e drithërave:

ema_data <- EMA(ts_prodhimi_dritherave, n = 3)

plot(ts_prodhimi_dritherave, main = "Seri Kohore dhe Exponential Moving Average", col = "blue", type = "l")
lines(ema_data, col = "purple", lwd = 2)
legend("topright", legend = c("Seri Kohore", "Exponential Moving Average"), col = c("blue", "purple"), lwd = 2)

EMA ka një peshë më të madhe për vlerat më të fundit të serisë kohore, duke reaguar më shpejt ndaj ndryshimeve sesa mesatarja e lëvizshme.
Vijat blu tregojnë serinë origjinale, ndërsa vijat vjollcë përfaqësojnë EMA.
EMA përshtatet më mirë me ndryshimet e fundit të prodhimit të drithërave. Për shembull, reagimi ndaj rritjes së shpejtë pas viteve 1990 është më i shpejtë në vijën vjollcë.
Është më e ndjeshme ndaj ndryshimeve të kohëve të fundit në krahasim me MA.
EMA është një mjet më i mirë kur është e rëndësishme të kapet sjellja e fundit e serisë kohore.
4.3 Sheshimi i thjeshtë eksponencial për prodhimin e drithërave:

data_ts <- ts(data$ProdhimiDritherave, start = c(1961), frequency = 1)

shesh_exp_drithera <- HoltWinters(ts_prodhimi_dritherave, beta = FALSE, gamma = FALSE)
shesh_exp_drithera

## Holt-Winters exponential smoothing without trend and without seasonal component.
## 
## Call:
## HoltWinters(x = ts_prodhimi_dritherave, beta = FALSE, gamma = FALSE)
## 
## Smoothing parameters:
##  alpha: 0.9999202
##  beta : FALSE
##  gamma: FALSE
## 
## Coefficients:
##       [,1]
## a 684021.6

plot(shesh_exp_drithera, 
     main = "Sheshimi i thjeshtë eksponencial për prodhimin e drithërave", 
     xlab = "Viti", 
     ylab = "Prodhimi i drithërave")

Ky grafik paraqet një model SES pa trend dhe sezonalitet.
Linjat e parashikimit tregojnë një mesatare të zbutur të serisë kohore.
Modeli SES përshtatet mirë me të dhënat, duke kapur trendin rritës në seri. Vlerat e alfa-s sugjerojnë një peshim të fortë të vlerave të fundit.
Nuk mund të trajtojë trendet ose sezonalitetet afatgjata, por është efektiv për seritë pa këto karakteristika.
SES është një metodë e thjeshtë për zbutjen e të dhënave dhe përshtatet mirë për të dhëna me variacion minimal.
4.4 Sheshimi eksponencial me Holt-Winters për prodhimin e drithërave:

# Aplikimi i modelit Holt-Winters me alpha dhe beta (pa gamma)
holt_model <- HoltWinters(data_ts, beta = 0.1, gamma = FALSE) # Gamma = FALSE heq sezonalitetin

# Shfaqja e parametrave të modelit
holt_model

## Holt-Winters exponential smoothing with trend and without seasonal component.
## 
## Call:
## HoltWinters(x = data_ts, beta = 0.1, gamma = FALSE)
## 
## Smoothing parameters:
##  alpha: 0.9999525
##  beta : 0.1
##  gamma: FALSE
## 
## Coefficients:
##         [,1]
## a 684022.202
## b   2815.481

# Vizualizimi i modelit Holt-Winters
plot(holt_model, 
     main = "Modeli Holt-Winters për Prodhimin e Drithërave (Me Trend)", 
     xlab = "Viti", 
     ylab = "Prodhimi i Drithërave")

Ky model përfshin trendin në mënyrë të qartë, por nuk përfshin sezonalitet. - Vijat e modeluara tregojnë një trend të lëmuar të rritjes së prodhimit të drithërave. - Modeli është shumë i përshtatshëm për të dhënat që tregojnë një trend rritës afatgjatë, si në këtë rast. - Parametrat alfa dhe beta përputhen mirë me të dhënat origjinale, duke rezultuar në një model të saktë. - Konkluzion: - Holt-Winters është një metodë më e avancuar për të trajtuar seri me trend dhe siguron një përshtatje më të mirë krahasuar me SES.

4.5 Modeli ARIMA për emëtimin e CO2:

library(forecast)
# Modeli ARIMA për Prodhimi i Drithërave
fit_eme_CO2 <- auto.arima(ts_CO2)
summary(fit_eme_CO2)

## Series: ts_CO2 
## ARIMA(1,0,0) with non-zero mean 
## 
## Coefficients:
##          ar1    mean
##       0.8695  1.5789
## s.e.  0.0583  0.2436
## 
## sigma^2 = 0.07664:  log likelihood = -7.76
## AIC=21.53   AICc=21.96   BIC=27.81
## 
## Training set error measures:
##                       ME      RMSE       MAE       MPE    MAPE    MASE
## Training set 0.006669262 0.2721802 0.1926131 -3.439636 13.6033 1.01612
##                     ACF1
## Training set -0.05593145

# Parashikimi me ARIMA
forecast_eme_CO2 <- forecast(fit_eme_CO2, h = 4)
autoplot(forecast_eme_CO2)

Ky grafik tregon modelin ARIMA(1,0,0) për emetimet e CO2, duke përfshirë vijën e parashikimit dhe kufijtë e besueshmërisë. - Modeli përputhet mirë me të dhënat origjinale të emetimeve të CO2. Autokorrelacionet e vonesave (lags) tregojnë se një urdhër i parë autoregresiv (AR) është i mjaftueshëm. - Kufijtë e besueshmërisë janë të ngushtë, duke sugjeruar një parashikim të besueshëm. - ARIMA është një metodë e fuqishme për të trajtuar seri kohore komplekse. Për CO2, modeli me komponent autoregresiv është i mjaftueshëm për kapjen e dinamikës së të dhënave.

4.6 Modeli ARIMA për sasinë e Plehrave Kimike:

library(forecast)
# Modeli ARIMA për Prodhimi i Drithërave
fit_plehra <- auto.arima(ts_plhera_kimike)
summary(fit_plehra)

## Series: ts_plhera_kimike 
## ARIMA(1,0,0) with non-zero mean 
## 
## Coefficients:
##          ar1      mean
##       0.9063  44410.14
## s.e.  0.0505  16696.52
## 
## sigma^2 = 196055624:  log likelihood = -657.8
## AIC=1321.6   AICc=1322.02   BIC=1327.88
## 
## Training set error measures:
##                    ME     RMSE      MAE       MPE     MAPE    MASE       ACF1
## Training set 938.7646 13766.64 8955.413 -19.40321 37.64379 1.07692 0.06175204

# Parashikimi me ARIMA
forecast_plehra <- forecast(fit_plehra, h = 4)
autoplot(forecast_plehra)

Ky grafik paraqet modelin ARIMA(1,0,0) për sasinë e plehrave kimike, së bashku me parashikimet për vitet e ardhshme. - Modeli kap rritjen e qëndrueshme të plehrave kimike gjatë periudhës. Autokorrelacionet sugjerojnë një trend të fortë autoregresiv për këtë variabël. - Vijat e parashikimeve tregojnë një trend rritës të vazhdueshëm, që është në përputhje me trendin historik. - ARIMA(1,0,0) është efektiv për trajtimin e kësaj serie, duke ofruar parashikime të besueshme për përdorimin e plehrave kimike.

5. Modelet Parashikuese

5.1 Ndarja e të dhënave në train dhe test:

data_train = head(data, 48) 
data_test = tail(data, 12)

5.2 Parashikimi me metodën Drift:

5.2.1) Ndërtimi i modelit dhe saktësia e tij

library(forecast)
library(fpp2)

## ── Attaching packages ────────────────────────────────────────────── fpp2 2.5 ──

## ✔ fma       2.5     ✔ expsmooth 2.3

##

drift.drithera <- rwf(data_train$ProdhimiDritherave, h = 12, drift = TRUE)
accuracy(drift.drithera, data_test$ProdhimiDritherave)

##                         ME     RMSE      MAE        MPE     MAPE      MASE
## Training set -1.888653e-11 95967.30 57358.58 -0.9262221 9.728369 0.9755229
## Test set      3.549491e+04 45834.02 38969.84  5.0902108 5.608772 0.6627773
##                   ACF1
## Training set 0.1057804
## Test set            NA

ME (Mean Error):
- Gabimi mesatar është shumë afër 0 në trajnim, që tregon një përshtatje të mirë të modelit me të dhënat historike.
- Në testim, ME rritet në 35,494.91, duke treguar një tendencë të vogël të modelit për të mbivlerësuar prodhimin e drithërave në të ardhmen.
RMSE (Root Mean Squared Error):
- Vlera RMSE në trajnim është 95,967.30, që sugjeron devijim të moderuar nga të dhënat aktuale. Në testim, RMSE është më i ulët (45,834.02), duke treguar parashikime relativisht të sakta.
MAPE (Mean Absolute Percentage Error):
- Në trajnim, MAPE është 9.73%, që tregon një gabim mesatar relativisht të ulët. Në testim, MAPE është 5.61%, duke treguar parashikime më të sakta për të ardhmen.

5.2.2) Kontrolli i mbetjeve të modelit

checkresiduals(drift.drithera)

## 
##  Ljung-Box test
## 
## data:  Residuals from Random walk with drift
## Q* = 8.081, df = 10, p-value = 0.6209
## 
## Model df: 0.   Total lags used: 10

Rezultatet nga grafiku i mbetjeve: - Grafiku i mbetjeve tregon se ato janë shpërndarë rastësisht rreth vijës zero, që është një tregues pozitiv se modeli nuk ka mbetje të strukturuara.

Autokorrelacioni i mbetjeve: - ACF për mbetjet tregon mungesë të autokorrelacioneve domethënëse, që sugjeron se mbetjet janë në mënyrë të pavarur nga njëra-tjetra dhe janë stacionare.
Testi Ljung-Box: - Rezultati i testit Ljung-Box tregon një p-vlerë prej 0.6209, që është mbi pragun 0.05. Kjo sugjeron se nuk ka autokorrelacion të rëndësishëm mes mbetjeve.

5.3 Parashikimi me metodën ARIMA:

5.3.1) Krijimi i modelit dhe rezultatet e tij

# Ngarkimi i bibliotekave
library(forecast)
library(ggpubr)

## 
## Attaching package: 'ggpubr'

## The following object is masked from 'package:forecast':
## 
##     gghistogram

# Krijimi i serisë kohore për trajnimin
ts_prodhimi_dritherave <- ts(data_train$ProdhimiDritherave, start = 1961, frequency = 1)

# Krijimi i modelit ARIMA me auto.arima
arima.drithera <- auto.arima(ts_prodhimi_dritherave)

# Rezultatet e modelit
summary(arima.drithera)

## Series: ts_prodhimi_dritherave 
## ARIMA(0,1,0) 
## 
## sigma^2 = 9.255e+09:  log likelihood = -605.98
## AIC=1213.95   AICc=1214.04   BIC=1215.8
## 
## Training set error measures:
##                    ME     RMSE      MAE       MPE     MAPE      MASE      ACF1
## Training set 6559.624 95193.06 57578.95 0.2119825 9.723723 0.9792708 0.1055468

Gabimet janë relativisht të ulëta (MAPE = 9.72%), që tregon një përshtatje të mirë të modelit me të dhënat historike.
RMSE (95,193.06) sugjeron një devijim të moderuar të modelit nga të dhënat aktuale.
Gabimet në testim janë më të larta, me një MAPE prej 11.41%, duke treguar një rënie të lehtë të saktësisë në parashikime për të ardhmen.
RMSE dhe MAE gjithashtu janë më të larta në testim, por priten për shkak të ndryshueshmërisë së të dhënave.

5.3.2)Parashikimi me modelin ARIMA dhe paraqitja grafike

# Parashikimi me modelin ARIMA
forecast_arima <- forecast(arima.drithera, h = 12)

# Vizualizimi i parashikimit
plot(forecast_arima, main = "Parashikimi me Modelin ARIMA", ylab = "Prodhimi i Drithërave", xlab = "Viti")
lines(ts(data_test$ProdhimiDritherave, start = 2009, frequency = 1), col = "red", lwd = 2)
lines(ts(arima.drithera$fitted, start = 1961, frequency = 1), col = "orange", lwd = 2)
legend("topleft", c("Train", "Fit", "Prediction", "Test"),
       fill = c("black", "orange", "blue", "red"), box.col = "black")

- Kufijtë e besueshmërisë janë relativisht të ngushtë, duke treguar besueshmëri të lartë në parashikime. - Të dhënat reale dhe parashikimet ndjekin një trend të ngjashëm, megjithëse ka disa devijime në vlerat specifike.

5.3.3) Saktësia e modelit ARIMA

accuracy(forecast_arima, data_test$ProdhimiDritherave)

##                     ME     RMSE      MAE        MPE      MAPE      MASE
## Training set  6559.624 95193.06 57578.95  0.2119825  9.723723 0.9792708
## Test set     78999.000 81525.74 78999.00 11.4104076 11.410408 1.3435711
##                   ACF1
## Training set 0.1055468
## Test set            NA

Performanca e modelit është më e mirë për setin e trajnimit se për setin e testimit, me disa metrika si RMSE, MAE, dhe MAPE që janë më të larta për setin e testimit, duke sugjeruar se modeli nuk është shumë i saktë në parashikimet për të dhënat që nuk ka parë më parë (setin e testimit).

5.3.4) Kontrolli i mbetjeve të modelit ARIMA

checkresiduals(arima.drithera)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(0,1,0)
## Q* = 8.1954, df = 10, p-value = 0.6098
## 
## Model df: 0.   Total lags used: 10

# Testi i Ljung-Box për autocorrelacionin e mbetjeve
Box.test(arima.drithera$residuals, lag = 10, type = "Ljung-Box")

## 
##  Box-Ljung test
## 
## data:  arima.drithera$residuals
## X-squared = 8.1954, df = 10, p-value = 0.6098

Veme re se autokorrelacionet e mbetjeve jane brenda intervalit te besimit dhe kjo do te thote se seria e mbetjeve eshte stacionare. Rezultatet sugjerojnë që nuk ka autokorrelacion të rëndësishme në mbetjet e modelit. P-vlera prej 0.6098 është më e madhe se 0.05, që nënkupton se mbetjet janë të pavarura dhe modelin ARIMA(0,1,0) mund të jetë i përshtatshëm për të dhënat.

5.3.5) Vizualizimi i shpërndarjes së mbetjeve

ggdensity(as.numeric(arima.drithera$residuals), fill = "lightgray", main = "Shpërndarja e Mbetjeve")

ggqqplot(as.numeric(arima.drithera$residuals), color = "magenta",main = "QQ-Plot i Mbetjeve")

Që modeli te jetë i mirë ne presim që kuantilet teorike dhe të zgjedhjes të jenë brenda zonës gri pra afer drejtëzës. Në këtë rast shumica e tyre janë brenda zonës gri.

5.4 Parashikimi me metodën RegARIMA me regresorë PlehraKimike:

5.4.1) Krijimi i modelit RegARIMA dhe rezultatet e tij

# Ngarkimi i bibliotekave të nevojshme
library(forecast)
library(ggpubr)
library(tseries)

# Krijimi i regresorit si seri kohore (PlehraKimike)
ts_plehra <- ts(data_train$PlehraKimike, start = 1961, frequency = 1)

# Krijimi i modelit RegARIMA me regresorë PlehraKimike
regarima_model <- auto.arima(
  ts_prodhimi_dritherave,
  xreg = cbind(ts_plehra),
  seasonal = FALSE
)

# Rezultatet e modelit
summary(regarima_model)

## Series: ts_prodhimi_dritherave 
## Regression with ARIMA(1,0,0) errors 
## 
## Coefficients:
##          ar1  intercept    xreg
##       0.7611  408293.30  4.7118
## s.e.  0.1053   50996.14  0.7259
## 
## sigma^2 = 5.179e+09:  log likelihood = -603.82
## AIC=1215.64   AICc=1216.57   BIC=1223.13
## 
## Training set error measures:
##                   ME     RMSE      MAE       MPE     MAPE      MASE       ACF1
## Training set 3379.92 69679.09 53657.94 -0.937631 9.133682 0.9125844 0.06813269

Ky është një model regresioni që përfshin gabime të ARIMA(1,0,0), dhe ka pasur performancë të mirë për setin e trajnimit, me metrika të sakta si ME, RMSE, dhe MAE që janë relativisht të vogla. AIC, AICc dhe BIC sugjerojnë se ky model mund të jetë një zgjedhje e mirë për këtë dataset.

5.4.2) Parashikimi dhe vizualizimi i tyre

# Regresorët për testimin
xreg_test <- cbind(data_test$PlehraKimike)

# Parashikimi për 12 vitet e ardhshme
forecast_regarima <- forecast(regarima_model, xreg = xreg_test, h = 12)

# Vizualizimi i parashikimeve
plot(forecast_regarima, main = "Parashikimi me Modelin RegARIMA",
     ylab = "Prodhimi i Drithërave", xlab = "Viti")
lines(ts(data_test$ProdhimiDritherave, start = 2009, frequency = 1), col = "red", lwd = 2)
lines(ts(regarima_model$fitted, start = 1961, frequency = 1), col = "orange", lwd = 2)
legend("topleft", c("Train", "Fit", "Prediction", "Test"),
       fill = c("black", "orange", "blue", "red"), box.col = "black")

Modeli parashikon mirë vlerat por me disa luhatje të vogla.

5.4.3) Saktësia e modelit RegARIMA

accuracy(forecast_regarima, data_test$ProdhimiDritherave)

##                    ME     RMSE      MAE       MPE     MAPE      MASE       ACF1
## Training set  3379.92 69679.09 53657.94 -0.937631 9.133682 0.9125844 0.06813269
## Test set     44722.91 56753.30 46993.34  6.432737 6.789821 0.7992366         NA

5.4.4) Kontrolli i mbetjeve të modelit

# Kontrollimi i rezidualeve
checkresiduals(regarima_model)

## 
##  Ljung-Box test
## 
## data:  Residuals from Regression with ARIMA(1,0,0) errors
## Q* = 9.2685, df = 9, p-value = 0.4129
## 
## Model df: 1.   Total lags used: 10

# Testi Ljung-Box për autocorrelacionin e mbetjeve
Box.test(regarima_model$residuals, lag = 10, type = "Ljung-Box")

## 
##  Box-Ljung test
## 
## data:  regarima_model$residuals
## X-squared = 9.2685, df = 10, p-value = 0.5068

Veme re se autokorrelacionet e mbetjeve jane brenda intervalit te besimit dhe kjo do te thote se seria e mbetjeve eshte stacionare. Rezultatet sugjerojnë që nuk ka autokorrelacion të rëndësishme në mbetjet e modelit. P-vlera prej 0.4129 është më e madhe se 0.05, që nënkupton se mbetjet janë të pavarura dhe modeli mund të jetë i përshtatshëm për të dhënat.

5.4.5) Shpërndarja e mbetjeve

# Shpërndarja e rezidualeve
ggdensity(as.numeric(regarima_model$residuals), fill = "lightblue", main = "Shpërndarja e Rezidualeve")

# QQ-Plot i rezidualeve
ggqqplot(as.numeric(regarima_model$residuals),color = "red" ,main = "QQ-Plot i Rezidualeve")

Që modeli të jetë i mirë ne presim që kuantilet teorike dhe të zgjedhjes të jenë brenda zonës gri pra afer drejtëzës. Në këtë rast shumica e tyre janë brenda zonës gri.

5.4.6) Paraqitja e gabimeve e Regresit dhe ARIMA

cbind("Gabimet e regresit" = residuals(regarima_model, type="regression"),
 "Gabimet ARIMA" = residuals(regarima_model, type="innovation")) %>%
 autoplot(facets=TRUE)

5.5 Parashikimi me një model Hibrid me regresorë të jashtëm PlehraKimike:

5.5.1) Ndërtimi i modelit

library(forecastHybrid)

## Loading required package: thief

hyb.mod.8 <- hybridModel(ts(data_train$ProdhimiDritherave,start=1961,frequency=1), models = "an",
 a.args = list(xreg = ts(data_train$PlehraKimike,start=1961,frequency=1)),
 n.args = list(xreg = ts(data_train$PlehraKimike,start=1961,frequency=1)))

## Fitting the auto.arima model

## Fitting the nnetar model

hyb.f8<-forecast(hyb.mod.8,xreg =ts(data_test$PlehraKimike,start=c(2009,1),frequency=1))

## Warning in forecast.hybridModel(hyb.mod.8, xreg = ts(data_test$PlehraKimike, :
## The number of rows in xreg should match h. Setting h to nrow(xreg).

autoplot(hyb.f8)+ylab("Prodhim Drithera")+
 autolayer(ts(data_test$ProdhimiDritherave,start=c(2009,1),frequency = 1),series="Data")

Modeli hibrid ndjek mirë trendin e prodhimit të drithërave gjatë periudhës së trajnimit. Vijat e parashikimit janë të përshtatshme dhe të njëtrajtshme me të dhënat reale. Parashikimet për periudhën e testimit janë të afërta me vlerat reale, duke treguar aftësi të mirë të modelit për të kapur dinamikën e të dhënave.

5.5.2) Rezultatet dhe saktësia e modelit

summary(hyb.mod.8)

##            Length Class          Mode     
## auto.arima 19     forecast_ARIMA list     
## nnetar     17     nnetar         list     
## weights     2     -none-         numeric  
## frequency   1     -none-         numeric  
## x          48     ts             numeric  
## xreg        2     -none-         list     
## models      2     -none-         character
## fitted     48     -none-         numeric  
## residuals  48     ts             numeric

accuracy(hyb.f8)

##                    ME     RMSE      MAE       MPE     MAPE      MASE       ACF1
## Training set 2763.768 57771.71 44621.36 -0.668358 7.425179 0.7588952 0.09308877

Vlerat e ulëta tregojnë një përshtatje të saktë të modelit me të dhënat historike. MAPE (7.38%): Gabimi mesatar relativisht i ulët tregon një model të saktë dhe të besueshëm për parashikimet. Përfshirja e plehrave kimike: Regresori i jashtëm përmirëson parashikimet duke ofruar informacion shtesë mbi ndikimin e jashtëm në prodhimin e drithërave

5.5.3) Kontrolli për mbetjet e modelit

checkresiduals(hyb.mod.8$residuals)

## 
##  Ljung-Box test
## 
## data:  Residuals
## Q* = 8.477, df = 10, p-value = 0.5824
## 
## Model df: 0.   Total lags used: 10

Shpërndarja e mbetjeve: - Mbetjet janë shpërndarë rastësisht rreth vijës zero, duke treguar se modeli nuk ka mbetje të strukturuara.

Autokorrelacioni i mbetjeve: - Grafiku ACF për mbetjet nuk tregon autokorrelacione domethënëse.

Testi Ljung-Box: - P-vlera e testit është 0.5903, që është mbi pragun 0.05, duke treguar mungesë të autokorrelacionit të mbetjeve.

5.5.4) Box Test për mbetjet e modelit

Box.test(x = hyb.mod.8$residuals)

## 
##  Box-Pierce test
## 
## data:  hyb.mod.8$residuals
## X-squared = 0.40728, df = 1, p-value = 0.5234

p-value = 0.5129 > 0.05, që do të thotë se nuk ka evidence për autokorrelacion në mbetje. Kjo sugjeron se modeli ka kapur informacionin e nevojshëm nga të dhënat dhe mbetjet janë të pavarura.

5.5.5) Paraqitja grafike e mbetjeve të modelit

library(ggpubr)
# Density plot
ggdensity(hyb.mod.8$residuals, fill = "lightgray")

## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_density()`).

# QQplot
ggqqplot(hyb.mod.8$residuals)

## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.
## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.
## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.

## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_qq()`).

## Warning: Removed 1 row containing non-finite outside the scale range (`stat_qq_line()`).
## Removed 1 row containing non-finite outside the scale range (`stat_qq_line()`).

Mbetjet nuk ndjekin një shpërndarje normale. Që modeli te jetë i mirë ne presim që kuantilet teorike dhe të zgjedhjes të jenë brenda zonës gri pra afer drejtëzës. Në këtë rast shumica e tyre janë brenda zonës gri.

5.5.6) Shapiro Test për mbetjet e modelit

shapiro.test(hyb.mod.8$residuals)

## 
##  Shapiro-Wilk normality test
## 
## data:  hyb.mod.8$residuals
## W = 0.95, p-value = 0.04332

Mbetjet nuk kanë shpërndarje normale.

5.6 Parashikimi me 4 modele të ndryshme:

5.6.1) Sheshimi i thjeshtë eksponencial dhe saktësia e tij

model1 <- ses(data_train$ProdhimiDritherave, h = 12)
accuracy(model1, data_test$ProdhimiDritherave)

##                     ME     RMSE      MAE        MPE      MAPE      MASE
## Training set  6553.061 95194.08 57574.23  0.2095708  9.722069 0.9791905
## Test set     79010.497 81536.88 79010.50 11.4120814 11.412081 1.3437667
##                   ACF1
## Training set 0.1056043
## Test set            NA

Modeli SES ka një MAPE prej 9.72% për setin e trajnimit dhe një 11.41% për setin e testit, duke treguar një përfundim të ngjashëm për të dyja periudhat. Ky është një model i thjeshtë, por mund të ketë vështirësi për të kapur dinamikën më komplekse të serisë kohore.

5.6.2) Sheshimi eksponencial me 2 parametra (metoda Holt) dhe saktësia e tij

# Sheshimi eksponencial me 2 parametra (metoda Holt)
model2 <- holt(data_train$ProdhimiDritherave, h = 12)
accuracy(model2, data_test$ProdhimiDritherave)

##                     ME     RMSE      MAE       MPE     MAPE      MASE
## Training set -13579.06 96733.31 55634.80 -3.006817 9.561125 0.9462058
## Test set      51189.36 56402.54 51189.36  7.370276 7.370276 0.8706002
##                    ACF1
## Training set 0.07601034
## Test set             NA

Modeli Holt ka një përformancë të mirë në setin e testimit (MAPE 7.37%) dhe ka përmirësuar performancën krahasuar me SES. Kjo tregon se modeli është në gjendje të kapë më mirë trendet në të dhëna.

5.6.3) Sheshimi eksponencial (me trend të zbutur) dhe saktësia e tij

# Sheshimi eksponencial (me trend të zbutur)
model3 <- holt(data_train$ProdhimiDritherave, damped = TRUE, h = 12)
accuracy(model3, data_test$ProdhimiDritherave)

##                     ME     RMSE      MAE       MPE      MAPE      MASE
## Training set -3198.034 93969.01 54998.30 -1.681557  9.496265 0.9353805
## Test set     73535.286 76158.92 73535.29 10.617319 10.617319 1.2506473
##                    ACF1
## Training set 0.08086838
## Test set             NA

Ky model ka një MAPE të lartë në setin e testimit (10.62%) dhe një RMSE të lartë, që sugjeron se modeli me trend të zbutur mund të jetë i dobishëm për disa lloje të serive kohore, por mund të kërkojë përshtatje të mëtejshme për të kapur më mirë tendencat.

5.6.4) Modeli ETS dhe saktësia e tij

# Modeli ETS
model4 <- forecast(ets(data_train$ProdhimiDritherave), h = 12)
accuracy(model4, data_test$ProdhimiDritherave)

##                      ME     RMSE       MAE        MPE     MAPE     MASE
## Training set   8326.474 100393.8  61572.48  0.1655844 10.53745 1.047191
## Test set     110880.276 112694.5 110880.28 16.0519018 16.05190 1.885790
##                   ACF1
## Training set 0.3126771
## Test set            NA

Modeli ETS ka një MAPE të lartë në setin e testimit (16.05%) dhe një RMSE më të madh në krahasim me Holt dhe SES. Kjo mund të sugjerojë se modeli është shumë kompleks dhe nuk është gjithmonë më i përshtatshëm për të gjitha llojet e serive kohore, duke pasur parasysh se mund të ndodhin të dhëna të paqëndrueshme

5.6.5) Vizualizimi i parashikimeve me këto modele

par(mfrow = c(2,2))
plot(model1)
plot(model2)
plot(model3)
plot(model4)

5.7 Parashikimi me modelin TBATS:

5.7.1) Ndërtimi i modelit

tbats.12<-tbats(ts(data_train$ProdhimiDritherave,start=1961,frequency=1))
tbats.f12<-forecast(tbats.12,12)
autoplot(tbats.f12)+ylab("Prodhimi i Dritherave")+
 autolayer(ts(data_test$ProdhimiDritherave,start=c(2009,1),frequency = 1),series="Data")

Modeli TBATS u përdor për të ndërtuar një model të serive kohore mbi prodhimin e drithërave. Modeli i ndërtuar u parashikua për 12 periudha (vite), dhe rezultatet u krahasuan me të dhënat reale të testit.Shikojmë një përshtatje të mirë të vijave.

5.7.2) Rezultatet dhe saktësia e modelit

summary(tbats.12)

##                   Length Class  Mode     
## lambda             1     -none- numeric  
## alpha              1     -none- numeric  
## beta               1     -none- numeric  
## damping.parameter  1     -none- numeric  
## gamma.values       0     -none- NULL     
## ar.coefficients    0     -none- NULL     
## ma.coefficients    0     -none- NULL     
## likelihood         1     -none- numeric  
## optim.return.code  1     -none- numeric  
## variance           1     -none- numeric  
## AIC                1     -none- numeric  
## parameters         2     -none- list     
## seed.states        2     -none- numeric  
## fitted.values     48     ts     numeric  
## errors            48     ts     numeric  
## x                 96     -none- numeric  
## seasonal.periods   0     -none- NULL     
## y                 48     ts     numeric  
## call               2     -none- call     
## series             1     -none- character
## method             1     -none- character

accuracy(forecast(tbats.12,12),data_test$ProdhimiDritherave)

##                     ME     RMSE      MAE       MPE      MAPE     MASE      ACF1
## Training set  5030.823 91887.58 60236.44 -1.300912 10.365632 1.024468 0.1295127
## Test set     21770.164 33000.51 28076.41  3.106058  4.044544 0.477508        NA

Përdorueshmëria e modelit për saktësi: MAPE = 4.04%, që tregon një gabim relativisht të vogël dhe model i përshtatshëm për këtë dataset.

5.7.3) Kontrolli i mbetjeve

checkresiduals(tbats.12$errors)

## 
##  Ljung-Box test
## 
## data:  Residuals
## Q* = 7.1862, df = 10, p-value = 0.7078
## 
## Model df: 0.   Total lags used: 10

Testi Ljung-Box: p-value = 0.7078, që tregon mungesë të korrelacionit të mbetjeve. Testi Shapiro-Wilk: p-value = 6.908e-05, duke treguar se mbetjet nuk ndjekin një shpërndarje perfekte normale. Veme re se autokorrelacionet e mbetjeve jane brenda intervalit te besimit dhe kjo do te thote se seria e mbetjeve eshte stacionare.

5.7.4) Box Test për mbetjet

Box.test(x = tbats.12$errors)

## 
##  Box-Pierce test
## 
## data:  tbats.12$errors
## X-squared = 0.16569, df = 1, p-value = 0.684

Një p-value më të madhe se 0.05 (në këtë rast, p = 0.684) që do të thotë se mbetjet janë të pavarura.

5.7.5) Vizualizimi i mbetjeve

library(ggpubr)
# Density plot
ggdensity(tbats.12$errors, fill = "lightgray")

## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.

# QQplot
ggqqplot(tbats.12$errors)

## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.
## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.
## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.

Vizualizimi përfshiu grafikun e densitetit dhe QQ-Plot, duke konfirmuar mungesën e shpërndarjes perfekte normale.Që modeli te jetë i mirë ne presim që kuantilet teorike dhe të zgjedhjes të jenë brenda zonës gri pra afer drejtëzës. Në këtë rast shumica e tyre janë brenda zonës gri.

5.7.6) Shaqpiro Test për mbetjet

shapiro.test(tbats.12$errors)

## 
##  Shapiro-Wilk normality test
## 
## data:  tbats.12$errors
## W = 0.86829, p-value = 6.908e-05

Një p-value shumë të vogël (në këtë rast, p = 6.908e-05) sugjeron që nuk mund të pranohet hipoteza e shpërndarjes normale.

7. Përfundime

Tendencat dhe Analiza e Serive Kohore: • Prodhimi i drithërave ka ndjekur një trend rritës në përgjithësi gjatë periudhës 1961-2020. • Variablat si përdorimi i plehrave kimike dhe emetimet e CO2 treguan ndikim të rëndësishëm në rritjen e prodhimit.
Korrelacioni midis Variablave: • Matrica e korrelacionit tregoi një lidhje të fortë pozitive midis plehrave kimike dhe prodhimit të drithërave (r = 0.89), si dhe midis CO2 dhe prodhimit të drithërave (r = 0.75). • Plehrat kimike rezultuan faktori kryesor që ndikon në prodhimin e drithërave.
Stacionariteti dhe Përgatitja e Modeleve: • Seria kohore për prodhimin e drithërave ishte jo-stacionare dhe u transformua përmes diferencimit për t’u bërë e përshtatshme për modelim. • Modelet u ndërtuan dhe krahasuan për të identifikuar qasjen më të përshtatshme për parashikime.
Modelet e Zbatuara: • Modeli i Regresit Linear tregoi ndikim të qartë të plehrave kimike në prodhimin e drithërave, me një R² = 0.79. • Modelet ARIMA, ETS, Holt-Winters dhe TBATS u përdorën për parashikime afatshkurtra dhe afatgjata. • Modeli TBATS tregoi performancë të mirë për shkak të fleksibilitetit të tij në trajtimin e trendit dhe sezonalitetit.
Performanca e Modeleve: • TBATS rezultoi më i saktë për dataset-in aktual me një gabim relativisht të ulët (MAPE = 4.04%). • Modelet si ETS dhe Holt shfaqën gabime më të larta krahasuar me TBATS dhe ARIMA.
Rezultatet Kryesore të Parashikimit: • Parashikimet treguan se prodhimi i drithërave do të vazhdojë të rritet, por rritja është e ndikuar ndjeshëm nga përdorimi i plehrave kimike. • Në rast të uljes së përdorimit të plehrave kimike ose emëtimeve të CO2, prodhimi mund të ketë ulje të konsiderueshme.
Vlerësimi i Gabimeve dhe Mbetjeve: • Analiza e mbetjeve tregoi që shumica e modeleve nuk kanë autokorrelacion të mbetjeve, duke treguar përshtatje të mirë me të dhënat. • Megjithatë, disa shpërndarje të mbetjeve devijuan nga normaliteti, duke sugjeruar nevojën për përmirësime në trajtimin e të dhënave ose përzgjedhjen e modeleve.

Detyra tregoi se modelet TBATS dhe ARIMA janë të përshtatshme për parashikimin e prodhimit të drithërave në bazë të serive kohore. Përmirësimi i vazhdueshëm i metodave dhe analizave mund të ndihmojë në marrjen e parashikimeve më të sakta për të ndihmuar politikat bujqësore dhe menaxhimin mjedisor.

8. Referencat

“Laboratore për Seritë Kohore” - Arnisa Sokoli
“Time Series Analysis and Its Applications: With R Examples” - Shumway & Stoffer
“Forecasting: Principles and Practice” - Hyndman & Athanasopoulos